#!/usr/bin/env perl
# by Tao Lang: langtao123xxx@126.com
# 
# 分类：
# 1. 以前的，还是单个单个比对；暂时还是先按以前的吧，跑起来再说。
# 2. 新的，只有一个 blastdb，使用 seqidlist
# 
# 2021/5/14 考虑了下，沿用以前的配置也行，然后全部使用 seqidlist。
# 
use strict;
use warnings;
use autodie qw/open close/;

use Getopt::Long;
use Cwd qw/getcwd abs_path/;
use File::Basename;
use JSON;

BEGIN {
    if (defined $ENV{MOJO_HOME}) {
        unshift @INC, "$ENV{MOJO_HOME}/lib";
    } else {
        unshift @INC, "./lib";
    }
}

use JsonConf qw/get_conf update_done load_done/;
use WriteData qw/write_data2file/;

# -------------------------------------------
# mod_perl
our $home_dir = $ENV{MOJO_HOME};
$home_dir = '.' unless (defined $home_dir and -d $home_dir);
# ===========================================
# some settings
our $conf_dir = join '/', $home_dir, 'conf/pp_locator';
our $conf_dir_species = "$conf_dir/species";
# 先读基础配置
our $conf = get_conf(join '/', $conf_dir, 'pp_locator.json');
our $data_dir = join '/', $home_dir, $conf->{data_dir};
our $public_tmp_dir = join '/', $home_dir, $conf->{public_tmp_dir};

# ===========================================
# command line args
# 
our $job_id;
my $help;

GetOptions(
    'jid|j=s' => \$job_id,
    'help|h' => \$help,
) or die "Error in command line arguments!\n";

check_opt();

# 收集信息
our $info = load_done($job_id, $public_tmp_dir);
# ===========================================
our $species = $info->{species};
load_conf_of_species($species);
our $chromosomes = $info->{chromosomes};
our $ofmt = exists $info->{ofmt} ? $info->{ofmt} : $conf->{ofmt};
my @specifiers = split /\s+/, $ofmt;
shift @specifiers;

# 染色体对应 blastdb 信息
our $bdb = bdb_info();

our $data_job_dir = "$data_dir/$job_id";
our $public_job_dir = "$public_tmp_dir/$job_id";
unless (-d $data_job_dir) {
    print STDERR "Error with 'jid': invalid!\n\n";
    exit 1;
}
our $fa_blastn = "$data_job_dir/blastn.fa";
our $fa_short = "$data_job_dir/blastn-short.fa";
# 中间文件
our $f_seqidlist = "$data_job_dir/seqid.list";

# -----------------------------------------------
# blastn
# 直接 F,R 比对到单条染色体上的各一个文件
my ($perc_identity, $qcov_hsp_perc) = @{$info}{qw/perc_identity qcov_hsp_perc/};
run_blastn();
$info->{status} = 'free'; # 运行结束，变更状态
update_done($info, $job_id, $public_tmp_dir);

# ===========================================
# subroutines
# ===========================================

sub check_opt {
    usage() if $help;

    unless (defined $job_id) {
        print STDERR "Show me your job id!\n\n";
        exit 1;
    }
}

sub usage {
    print STDERR <<USAGE;
For pp_locator web service.

Usage:
    $0 -j job-id

Options:

    -j, --jid
        Job id, like 't607440fb'.
        Required.

    -h, --help
        Display this help and quit.
        Optional.

USAGE
    exit 0;
}

# 比如对于短柄草 1 => 1Bd，使得更通用
sub bdb_info {
    our($species,$conf);
    my %h;

    my $href = $conf->{species}{$species}{blastdb}{chromosome};
    for my $g (keys %$href) {
        for my $n (keys %{ $href->{$g} }) {
            @{ $h{ $n . $g }}{qw/dbname sseqid/} = @{ $href->{$g}{$n} };
        }
    }

    return \%h;
}

# 将多个序列逐次 blastn 至各个染色体 blastdb
# 不返回啥

# 搞 2 个文件，1 个普通，一个 short
# 多个 sid 一起比对一个 chr
# 
# 直接处理后输出到统一的 blast 结果文件，使用 1A, 2A 这样的 sseqid
# 
sub run_blastn {
    my $db_dir = $conf->{species}{$species}{blastdb}{db_dir};
    my $n_threads = $conf->{threads};
    
    my $header_line = join("\t", @specifiers) . "\n";
    # 需要确定 blastdb 真实的染色体 id
    for my $cname (@$chromosomes) {
        my $f_out = "$public_job_dir/$cname.tsv";
        # write_data2file(\$header_line, $f_out);
        my $db = join("/", $db_dir, $bdb->{$cname}{dbname});
        my $sseqid = $bdb->{$cname}{sseqid};
        write_data2file(\$sseqid, $f_seqidlist);
        my $command;
        if (-T $fa_blastn) {
            my ($evalue, $word_size) = @{$info}{qw/evalue word_size/};
            $command = "blastn -num_threads $n_threads -query $fa_blastn -task blastn -db $db -outfmt '$ofmt' -perc_identity $perc_identity -qcov_hsp_perc $qcov_hsp_perc -evalue $evalue -word_size $word_size";
        }
        if (-T $fa_short) {
            my ($evalue_short, $word_size_short) = @{$info}{qw/evalue_short word_size_short/};
            $command = "blastn -num_threads $n_threads -query $fa_short -task blastn-short -db $db -outfmt '$ofmt' -perc_identity $perc_identity -qcov_hsp_perc $qcov_hsp_perc -evalue $evalue_short -word_size $word_size_short";
        }
        $command .= " -seqidlist $f_seqidlist";
        $command .= " -dust no" unless ($info->{dust});

        open my $fh_i, "$command |";
        open my $fh_o, '>', $f_out;
        print $fh_o $header_line;
        while (<$fh_i>) {
            next if /^\s*$/;
            my %tmp;
            @tmp{@specifiers} = split /\t/;
            $tmp{sseqid} = $cname;
            print $fh_o join("\t", @tmp{@specifiers});
        }
        close $fh_i;
        close $fh_o;
    }
}

sub load_conf_of_species {
    my $species = shift;
    my $f = join '/', $conf_dir_species, "$species.json";
    $conf->{species}{$species} = get_conf($f);

    return $conf;
}